addpath(genpath(pwd))
clear
% **************************** USER SETTING *******************************
orderApp = 3;
unitFree = 1;

%% Model specific part
[f,x,xp,y,yp,symparams,Phi] = DSGEmodel(unitFree);

% A string with the name of the steady state file
nameSteadyStateFile = '[auxOut,errorMes] = DSGEmodel_ss(params);';

%% The analytical derivatives of the model
% Folder were files are stored
positionFiles     = [pwd,'\ModelSpecification'];

% The parameters of the model
modelParams = cell(1,length(symparams));
for i=1:length(symparams)
    modelParams{i} = char(symparams(i));
end

% Function to evaluate the f-function in the steady state
strucOfArrays = struct('f',f);
nameOfFunction = 'numF.m';
generateMfuncNum_file(nameSteadyStateFile,modelParams,y,x,xp,yp,strucOfArrays,positionFiles,nameOfFunction)

% Derivatives of the model
setupModel=differentiate_dsge(f(1:end-length(Phi),1),yp,y,xp,x,symparams,orderApp,Phi);
setupModel.symparams = symparams;
save([positionFiles,'\setupModel'],'setupModel') % the structure model will be needed when we solve the model.
    
%% For the nominal bond prices
nameOfFunction = 'numDerivM.m';
[M,Mx,Mxp,My,Myp,Mxx,Mxxp,Mxy,Mxyp,Mxpx,Mxpxp,Mxpy,Mxpyp,Myx,Myxp,Myy,Myyp,Mypx,Mypxp,Mypy,Mypyp] ...
    = Anal_PricingKernel_derivatives(x,xp,y,yp,orderApp,1);
strucOfArraysMkernel = struct('M',M,'Mxp',Mxp,'Myp',Myp,...
    'Mxpx',Mxpx,'Mxpxp',Mxpxp,'Mxpy',Mxpy,'Mxpyp',Mxpyp,...
    'Mypx',Mypx,'Mypxp',Mypxp,'Mypy',Mypy,'Mypyp', Mypyp);
generateMfuncMkernel3rd(nameSteadyStateFile,modelParams,y,x,xp,yp,strucOfArraysMkernel,positionFiles,nameOfFunction)

%% For the real bond prices
nameOfFunction = 'numDerivM_Real.m';
[M,Mx,Mxp,My,Myp,Mxx,Mxxp,Mxy,Mxyp,Mxpx,Mxpxp,Mxpy,Mxpyp,Myx,Myxp,Myy,Myyp,Mypx,Mypxp,Mypy,Mypyp] ...
    = Anal_PricingKernel_Real_derivatives(x,xp,y,yp,orderApp,1);
strucOfArraysMkernel = struct('M',M,'Mxp',Mxp,'Myp',Myp,...
    'Mxpx',Mxpx,'Mxpxp',Mxpxp,'Mxpy',Mxpy,'Mxpyp',Mxpyp,...
    'Mypx',Mypx,'Mypxp',Mypxp,'Mypy',Mypy,'Mypyp', Mypyp);
generateMfuncMkernel3rd(nameSteadyStateFile,modelParams,y,x,xp,yp,strucOfArraysMkernel,positionFiles,nameOfFunction)
